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ABSTRACT 

We present high-cadence observations and simulations of the solar photosphere, obtained using the 
Rapid Oscillations in the Solar Atmosphere imaging system and the MuRAM magneto-hydrodynamic 
code, respectively. Each dataset demonstrates a wealth of magneto-acoustic oscillatory behaviour, 
visible as periodic intensity fluctuations with periods in the range 110 - 600 s. Almost no propagating 
waves with periods less than 140 s and 110 s are detected in the observational and simulated datasets, 
respectively. High concentrations of power are found in highly magnetised regions, such as magnetic 
bright points and inter granular lanes. Radiative diagnostics of the photospheric simulations replicate 
our observational results, confirming that the current breed of magneto-hydrodynamic simulations 
are able to accurately represent the lower solar atmosphere. All observed oscillations are generated 
as a result of naturally occurring magnetoconvective processes, with no specific input driver present. 
Using contribution functions extracted from our numerical simulations, we estimate minimum G- 
band and 4170 A continuum formation heights of 100 km and 25 km, respectively. Detected magneto- 
acoustic oscillations exhibit a dominant phase delay of —8° between the G-band and 4170 A continuum 
observations, suggesting the presence of upwardly propagating waves. More than 73% of MBPs (73% 
from observations, 96% from simulations) display upwardly propagating wave phenomena, suggesting 
the abundant nature of oscillatory behaviour detected higher in the solar atmosphere may be traced 
back to magnetoconvective processes occurring in the upper layers of the Sun's convection zone. 

Subject headings: methods: numerical — magnetohydrodynamics (MHD) — Sun: atmosphere — Sun: 
oscillations — Sun: photosphere 



1. INTRODUCTION 

Following the discovery of solar oscillations in the 1960s 
(|Leightonlll96Q[ h there has been a multitude of observa- 
tional evidence brought forward verify ing the existence of 
wave motion in the solar atmosphere (|Stein fc Leibacherl 
[1971 . In more recent years, the role of these oscillations 
is often closely linked to the search for efficient heat- 
ing mechanisms of the upper solar atmosphere (see e.g. 

traus et al.1l2008l : iCauzzi et al"1l2009l : iTarovan fe Erdelvil 
200C 



Wave phenomena in the lower solar atmosphere (pho- 
tosphere and chromosphere) can exist in a variety of 
different structures, and importantly exhi bit various 
modes of oscillation (jEdwin fc Roberts! 1 1983[ ). Magneto- 
acoustic, Alfven, sausage, and kink modes have re- 
cently been observed in a wide range of solar fea- 
tures dDe Poritieu et al ] 120071 : Uess et al.l l2007al 120091 : 
iMorton et al.ll2011l ). Importantly, however, is the fact 
that the majority of these wave phenomena exist in 
structures which are highly magnetic in nature, such 
as active regions, pores, and magnetic bright points 
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(MBPs). Indeed, iRimmele et~aT1 (|1995f ) found photo- 
spheric oscillations to preferentially exist in intrinsically 
magnetic intergranular lanes. In the Sun's atmosphere, 
magnetic field lines clump into tight bundles forming 
highly- structured wave guides, which allow oscillatory 
motion to be readily channelled along t hese structures 
(jShelvag et al.l l2010l : iFedun et al.l l2011al ). At any one 
time it is estimated that up to 2.2% of the entire so- 
lar su rface is covered by MBPs (Sanc hez Almeida et al.l 
120101 ). and with their intr insic magnetic fie ld strengths 
exceeding one kiloGauss (jJess et al.l l2010al h they may 
provide the basis of a suitable mechanism for efficient 
energy transport into the higher layers of the Sun's at- 
mosphere. 

Acoustic events occurring in the lower solar atmo- 
sphere are often closely examined in an attempt to de- 
termine the mechanisms behind, and t he energy fluxes 
associ ate with, such wave phenomena. Hoekz ema et al.l 
(2002) utilised high-resolution ground-based observa- 
tions to investigate the relation between upwardly 
propagating acoustic wa ves and co-spatial chr omo- 
spheric brightenings, while IStein fc Nordlundl (jl991[ ) sug- 
gested how acoustic modes may be excited stochas- 
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Fig. 1. — Simultaneous ROSA images acquired through G-band (left) and 4170 A continuum (right) niters at 14:06:21 UT on 2009 May 
28. Axes are in heliocentric arcseconds. Note the abundance of MBPs, particularly visible as intensity enhancements in the G-band. 



tically by non-adiabatic fluctuations in the gas pres- 
sure near the solar surface. Hi ghly energetic acous- 
tic wav es have been detected by iBello Gonzalez et al.l 
(|2010bl ) using spectropolarimetric data taken with the 
IMaX/ SUNRISE balloon experiment, with energy fluxes 
exc eeding 7000 Wm~ 2 i n the low photosphere. Contrar - 
ily, iFossum fc Carlssonl (|2005ft and lCarlsson et al.l (|2007ft 
examined the presence of acoustic oscillations in the so- 
lar chromosphere, and concluded that the energy flux 
of these waves was not sufficiently high to constitute 
the dominant atmospheric heating mechanism. How- 
ever, due to telemetry restrictions, the use of data ob- 
tained from space-based instrumentation limited their 
search to oscillations with periods longer than 25 s. This 
may have resulted in a considerable fractio n of energy, at 
short periods (~10 s), being ove rlooked (jDeForestl 120041 : 
lHasan fc van Ballegooijenll2008ft . Moreover, the contri- 
bution f unctions associated with t he fil ter bandpasses 
used by IFossum fc Carlssonl (|2005l ) and iCarlsson et al.l 
(|2007h are very extended, and must therefore be taken 
into further consideration (jBello Gonzalez et al.l 12009). 
Thus, the direct role of acoustic waves i n the heating of 
the solar a tmosp here remains uncerta in jK alkofen 2007; 
iBeck et alJfeOQa IBello Gonzalez et al.ll2010aft . To under- 
stand the presence of wave phenomena in the Sun's outer 
atmosphere, it is imperative to study the photospheric 
counterpart of these oscillations, and establish how ubiq- 
uitous these fundamental waves are in the lower solar at- 
mosphere. In this paper, we utilise a high cadence multi- 
wavelength data set to investigate propagating wave phe- 
nomena detected in observations and simulations of the 
lower solar atmosphere. 

2. DATA SETS 

2.1. Observations 

The observational data presented here are part of a se- 
quence obtained during 13:46 - 14:47 UT on 2009 May 
28, with the Dunn Solar Telescope at Sacramento Peak, 
New Mexico. We employ the Rapid Oscil lations in the 
Solar Atmosphere (ROSA; iJess et aT1l2010bf ) six-camera 



system to image a 69.3" x 69.1" region positioned at solar 
disk center. A spatial sampling of 0.069" per pixel was 
used for the ROSA cameras, to match the telescope's 
diffraction-limited resolution in the blue continuum to 
that of the CCD. This results in images obtained at 
longer wavelengths being slightly over sampled. However, 
this was deemed desirable to keep the dimensions of the 
field-of-view the same for all ROSA cameras. 

During the o bservations, high-order adaptive optics 
(jRimmeld 12004ft were used to correct wavefront defor- 
mations in real-time. The acquired images were further 
improved throug h use of the KISIP speckle reconstruc- 
tion algorithm (Wo ger et al.l 12008ft . utilising 16 — » 1 
restorations. Combining consistently excellent atmo- 
spheric seeing conditions with speckle image reconstruc- 
tion, a 2-pixel diffraction-limited resolution of ^100 km 
is achieved for both G-band and 4170 A blue continuum 
imaging. For the purposes of this paper, only obser- 
vations acquired through the photospheric G-band and 
4170 A blue continuum filters will be presented, each 
with a post-reconstruction cadence of 0.528 s. The G- 
band and 4170 A blue o continuum filters have central 
wavelengths of 4305.5 A and 4170.0 A, in addition to 
full-width half-maximums of 9.2 A and 52.0 A, respec- 
tively. To insure accurate coalignment in all bandpasses, 
broadband time series were Fourier co-registered and de- 
stretched using a 40 x 40 grid, e quating to a ~l-7" sepa - 
ration between spatial samples (|Jess et al.ll2007ai 12008). 

2.2. Simulations 

We utilise the MuRAM (jVogler et al.l 12005ft radia- 
tive magneto-hydrodynamic (MHD) code to investigate 
whether well-understood magnetoconvection processes 
can produce comparable phenomena to those detected 
in our observational time series. This code solves large- 
eddy radiative three-dimensional MHD equations on a 
Cartesian grid, and employs a fourth-order Runge-Kutta 
scheme to advance the numerical solution in time. The 
numerical domain has a physical size of 12 x 12 Mm 2 in 
the horizontal direction, 1.4 Mm in the vertical direction, 
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Fig. 2. — Upper Left: Histograms of the number of pixels, as a percentage of the total number of pixels, versus their normalised intensity 
in the G-band (black) and 4170 A continuum (blue) nltergrams, evaluated over the entire 61 minute duration of the data set (7.0xl0 9 
individual pixels). Upper Right: A zoom in of the intensity occurrences contained within the red box in the left panel. G-band intensities 
display an extended tail, well above the 4170 A continuum maximum, as a result of the presence of MBPs. Lower Left: Occurrence of 
intensities for isolated MBP structures, clearly displaying excess high-intensity measurements for G-band observations, also apparent in 
the upper-right panel. Lower Right: Intensity histograms for regions of the field-of-view not containing MBP structures. A narrowing of 
the 4170 A continuum profile near its maximum occurrence, caused by the removal of MBP intensities, directly results in the reduction of 
its positive transmission skewness. 



and is resolved by 480 x 480 x 100 grid cells, respectively. 
Our starting point for the simulations is a well-developed 
non-magnetic (B = 0) snapshot of photospheric con- 
vection taken approximately 2000 s (about 8 convec- 
tive turnover timescales) from the initial plane-parallel 
model. A uniform vertical magnetic field of 200 G was 
introduced at this stage, and a sequence of 872 snapshots 
recorded, each separated by a time interval of ~2 s. The 
resulting sequence, which is used for further radiative di- 
agnostics, covers approximately 30 min of physical time, 
corresponding to ~4 - 6 granular lifetimes. Full details 
of the simulatio n parameters and dom ain characteristics 
can be found in lShelvag et al.l (|2011bl ). 

The broadband filters used in the optical observations 
allow us to neglect the effects of line-of-sight velocities 
and magnetic fields in the radiative diagnostics. Due to 
a lack of detailed information on absorption line profiles 
contained within the filter bandpasses, precise compu- 
tation of the wavelength dependence of emergent inten- 



sity within each filter is not necessary. Thus a simplified 
and computationally efficient approach of calculating the 
transmitted filter inten sities is used. The methodology, 
similar to that used by Afran Tet al.l ([20111), is based o n 
the ATL AS9 spectral synthesis package of lKuruc 
Continuum and line opacity tables were calculated for a 
realistic solar chemical composition, incorporating a mi- 
croturbulent velocity, v tU rb — 2 kms -1 , and subsequently 
extracted from the ATLAS9 package using the subset 
codes XNFDF and DFSYNTHE. This procedure takes 
into account absorption by diatomic CH molecules. To 
avoid some of the problems highlighted by lAfram et al.l 
(2011), the line opacity wavelength grid between 4290 - 
4400 A was made denser, and centered at the wavelength 
of the G-band filter tran smission maximum. The DELO 
solver (|Rees et al.l 119891 ) was used to numerically inte- 
grate the corresponding radiative transport equations for 
all 480 x 480 vertical rays of the 872 photospheric snap- 
shots, allowing the computation of the intensity spectra 



4 



D.B. Jess et al. 






Fig. 3. — A G-band intensity image, averaged over the entire 61 minute duration of the time series (left). Green contours (middle 
and right panels) highlight the locations where waves are observed to propagate (i.e. a non-zero phase shift) between the G-band and 
4170 A continuum bandpasses. These contours are generated irrespective of oscillation frequency, and so cannot be used to pinpoint where 
particular periodicities exist. Dashed white lines outline a sub-field, which is magnified in the right panel. Blue and red contours in the 
right panel represent propagating wave power which is four and six orders-of-magnitude, respectively, above the background. A preferential 
existence of propagating waves is found in locations of increased magnetic field strength, such as MBPs. 



within the 4170 A and G-band filter bandpasses. Cal- 
culated spectra were then convolved with the G-band 
and 4170 A filter transmission functions to obtain a true 
intensity time series, which is directly comparable with 
our observational data. A more in-depth dis cussion of the 
method ology outlined here can be found in Shel vag et al.l 
([20llch . 

A natural feature of the numerical simulation is the 
ability to determine wavelength-dependent contribution 
functions of emergent radiation for the filter bandpasses 
which reside within its 1.4 Mm vertical domain. Us- 
ing previously obtained continuum intensities, as well as 
continuum and line band absorption coefficients, we de- 
termined the corresponding continuum and "band de- 
pression" (defined as (I — I c ) /I c , where / and I c are the 
emergent and continuum intensities, respectively) contri- 
bution functions for our observational G-band (4305.5 A 
central wavelength, 9.2 A bandpass) and 4170 A contin- 
uum (4170.0 A central wavelength, 52.0 A bandpass) fil- 
ters. The band depression contribution functions were 
calculated from the Rosseland means of the continuum 
and line absorptio n coefficients, and based on the defi- 
nition provided by iGurtovenko et al.l ([19741 ) and iMagainl 
(CHi). 

3. ANALYSIS AND DISCUSSION 

Overall, the Sun was very quiet on 2009 May 28, with 
no active regions, pores or large-scale magnetic activity 
present on disk. However, our ROSA field-of-view, posi- 
tioned at the center of the solar disk, contained numerous 
collections of MBPs, providing an ideal opportunity to 
examine these kiloGauss structures without the line-of- 
sight effects associated with off disc-center observations. 
Figure [1] shows simultaneous G-band and 4170 A contin- 
uum snapshots, revealing a wealth of MBPs, which are 
particularly visible as intensity enhancements in the G- 
band filtergram. 

3.1. Filter Transmission Profiles 

Following speckle image reconstruction, and to assist 
identification of inherently bright and dim structures, 
each of the G-band and 4170 A continuum images were 
normalised to their mean. To examine the throughput 



of each filter, an intensity histogram was created by bin- 
ning pixel intensities into intervals of 0.05 normalised 
units. The upper panels of Figure [2] plot histograms of 
the number of pixels (as a percentage of the total number 
of pixels) versus normalised intensity, summed over the 
entire 61 minute duration of the data set (in excess of 
7.0 x 10 9 individual pixels). It is clear from the upper- 
left panel of Figure [2] that the G-band mean is close to 
its modal value. However, for the 4170 A continuum ob- 
servations, the modal value of the images lies below the 
filter average, suggesting a degree of positive skew in the 
transmitted intensities. A measurement of the skewness 
can be determined through the calculation of Fisher and 
Pearson coefficients of skewness, given by, 



Fisher Coefficient = 



3/2 



Pearson Coefficient 



(^Ef=i(^-s) 2 ) 

3(*E£i(*«) 



(^£f = i(^- 50 2 ) 



and (1) 



(2) 



1/2 



where X{ is the pixel intensity, x is the mean, and x is the 
median, calculated over N independent measurements. 
While a symmetric distribution will yield a skewness co- 
efficient equal to zero, Fisher and Pearson statistics for 
the 4170 A continuum image sequence yield skewness co- 
efficients of 0.39 and 0.36, respectively. A potential weak- 
ness of the Fisher coefficient is its sensitivity to anoma- 
lous measurements occurring at the extreme ends of the 
intensity scale. However, consistency with the more ro- 
bust Pearson coefficient indicates the presence of a well 
defined, and positively skewed, 4170 A continuum trans- 
mission profile. 

Through examination of Figure [TJ it is clear that MBPs 
do not exhibit the same large contrast in 4170 A con- 
tinuum images that they do in simultaneous G-band 
snapshots. This can also be verified from the upper- 
right panel of Figure [2j where the occurrence of 4170 A 
continuum intensities falls to zero beyond 1.7 times the 
mean, compared with G-band intensities often exceeding 
twice the image mean. Interestingly, an intensity his- 
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Fig. 4. — Occurrence of oscillations simultaneously visible in both G-band and 4170 A continuum images, as a function of oscillation period 
and phase angle for MBP (left), and non-MBP (right) regions. The upper panels represent oscillations occurring in the observational time 
series, while the middle and lower panels display oscillatory phenomena detected in the noise- free and 5% added noise simulated datasets, 
respectively. The color scale represents the number of detections as a percentage of the total events with an associated coherence level 
exceeding 0.85. Consequently, the sum of each panel equals 100%. Dotted lines highlight a region inside which detections become unreliable 
due to cadence restrictions (0.528 s and 2 s for observational and simulated time series, respectively), while a horizontal dashed line represents 
a phase angle of degrees. A negative phase shift implies that an oscillation detected in the 4170 A continuum leads one also detected in 
the G-band. 



togram of isolated MBP structures (lower-left panel of 
Fig. [2]) reveals that the most commonly occurring MBP 
intensity is identical («1.1 times the image mean) in 
both bandpasses. However, an extended high- intensity 
tail is pronounced in the G-band MBP histogram, re- 
sulting in a high feature contrast for these structures. 
This elevated contrast is due to its continuum forma- 
tion level being depressed into deeper (and hotter) layers 
of the photosphere, where stronger magnetic fields and 
reduced CH abundances exist. An extended upper in- 
tensity range, and subsequently high feature contrast, of 
G-band MBP observations makes them ideal for the de- 
tection and tracking of such highly magnetic structures. 
As a result, we applied the MBP d etection an d tracking 
algori thms of lCrockett et al.l (|2QQ9Ll2Q10f ) and lKevs etaLl 
(|2Q1 lh to create binary maps for each of the G-band im- 
ages. These binary images are used to pinpoint the ex- 
act locations within our field-of-view where MBPs exist, 



thus allowing their intensity fluctuations to be precisely 
examined as a function of time. 

To investigate whether a shift of MBP intensities to 
lower values directly contributes to the positive skewness 
observed in the 4170 A continuum transmission profile, 
the MBP binary maps were utilised to mask out loca- 
tions containing these small-scale structures. The lower- 
right panel of Figure [2] displays the resulting intensity 
histogram for non-MBP areas. It is clear that the ex- 
tended G-band high-intensity tail is no longer present, 
in addition to a substantially reduced profile width near 
the 4170 A continuum maximum. Determination of the 
Fisher and Pearson statistics for the 4170 A continuum 
image sequence containing no MBPs yield skewness co- 
efficients of 0.20 and 0.18, respectively. These values are 
significantly reduced from those derived from the entire 
field-of-view, indicating that the lower MBP intensities 
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Fig. 5. — The normalised filter contributions as a function of atmospheric height for G-band (solid line) and 4170 A continuum (dashed line) 
bandpasses. The upper-left panel displays the contribution functions for a typical MBP, with the upper-right panel presenting contributions 
arising from the center of a granular cell. The bottom panel represents the contribution functions following a spatial averaging of the entire 
simulated field-of-view. A peak separation is clearly visible in the MBP case, while a peak overlap is present for the granular cell. However, 
due to an extended tail in all G-band contribution functions, a measurement of the center-of-gravity separation still yields a height difference 
in the case of a granular cell, even when the contribution peaks overlap. The atmospheric heights corresponding to the center-of-gravity of 
the G-band and 4170 A continuum contribution functions, in addition to the location of the visible solar surface, are represented in each 
panel by solid, dashed and dotted vertical lines, respectively. A 2-dimensional representation of these methods are shown in Figure [6] 



found in 4170 A continuum observations directly con- 
tributes to the positive skew found in the transmission 
profile. 

3.2. Propagating Wave Phenomena 

Following the methodology of Uess et al. ( 2007a. b]), the 
wavelet analysis routines of iTorrence fc Compol (j!998) 
were applied to both the G-band and 4170 A continuum 
time series to search for the presence of oscillatory be- 
haviour. Considerable oscillatory power, in both band- 
passes, was present throughout the entire field-of-view, 
with enhanced concentrations at the locations of MBPs. 
Due to the image plate scales and camera dimensions 
being identical for G-band and 4170 A continuum data, 
no further co-alignment or de-stretching between data 
sets was required to provide accurate structural align- 
ment. This can be verified visually through the exam- 
ination of Figure [TJ Subsequently, to test the propaga- 
tion of oscillatory phenomena between the two filters, 
coherence and phase information between simultaneous 
G-band and 4170 A continuum pixels were extracted. 
The resulting data array contained oscillation period 
(1.056 {Nyquist} — > 1830 seconds), power as a percent- 
age above the background, coherence (0.0 for incoherent 
waves — >> 1.0 for coherent oscillations), and phase angle 



(-180° +180°) for each of the 7.0xl0 9 imaging pixels. 
Following this nomenclature, a negative phase angle can 
be interpreted as an oscillation observed in the 4170 A 
continuum leading one also detected in the G-band. A 
time- averaged G-band image is shown in the left panel of 
Figure [3j contoured, in the middle panel, with the loca- 
tions of propagating (i.e. a non-zero phase shift) waves. 
Approximately 1% of the observational field-of-view con- 
tains propagating waves, with a preferential existence of 
these oscillations in locations of increased magnetic field 
strength, such as MBPs (right panel of Fig. [3j). 

To evaluate whether a preferential phase angle existed 
between G-band and 4170 A continuum oscillations, an 
occurrence plot, as a function of period and phase angle, 
was generated for all periodicities which displayed coher- 
ence values above 0.85. To differentiate between oscilla- 
tions occurring in MBP structures from those occurring 
in the intergranular lanes, the MBP binary maps created 
above were utilised. By multiplying through each period, 
power, coherence, and phase angle value by the corre- 
sponding pixel binary, an occurrence plot based on where 
the oscillation was generated can be displayed. The up- 
per panels of Figure |4] show the occurrence distributions 
for oscillations generated in both MBP and non-MBP 
regions. It is clear that a negative phase shift exists in 
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both panels, with peak occurrences at a phase angle of 
—8°. Positive phase angles do exist, suggesting G-band 
oscillations lead their 4170 A continuum counterparts in 
these instances. However, a clear majority of detections 
(73% for MBP structures and 71% for non-MBP regions) 
exhibit negative phase shifts, implying a preference for 
4170 A continuum oscillations to lead their G-band coun- 
terparts. A non-zero phase shift between the detected 
oscillations suggests that the observed wave phenomena 
are generated at distinctly different heights in the so- 
lar atmosphere, and implies the presence of propagating 
waves. Previous observational and theoretical studies of 
propagating wave phenomena have indicated their poten- 
tial for transporting energy away from the solar surface, 
and into higher atmosphe ric layers (e.g. Carlsso n fc Steinl 
[19971 : IVecchio et aLll2QQ7h . 

The simulated G-band and 4170 A continuum time se- 
ries, first interpolated on to a constant-cadence grid, 
were subjected to the same rigorous wavelet and MBP 
detection routines as described above. Following the 
evaluation of period, power, coherence, and phase an- 
gle between the simulated data, we are able to create a 
simulation-specific occurrence plot as a function of pe- 
riod and phase angle, as for the observational data. The 
middle panels of Figure |4] show the occurrence distribu- 
tions for oscillations generated in simulated MBP (left) 
and non-MBP (right) regions. As with the observational 
time series (upper panels of Fig. |4]), it is clear that a 
negative phase shift appears preferential. The peak oc- 
currences are consistent with the observational data, at 
a phase angle of —9°. In the simulated occurrence plot, 
96% of MBP structures display a negative phase shift, 
compared with 95% for non-MBP regions. Also of in- 
terest is the number of pixels undergoing oscillations in 
both the observational and simulated time series. Incor- 
porating both MBP and non-MBP oscillations, a total 
of 6.06 x 10 7 observational pixels (0.8% of the 7 x 10 9 
total pixels) are observed to oscillate, compared with 
1.62 x 10 7 simulated pixels (8.1% of the 2 x 10 8 total 
pixels). The higher occurrence rate of oscillations in the 
simulated data is most likely a result of the noise-free 
time series. Detector dark current, readout noise, and 
processing abnormalities will reduce the coherence levels 
of oscillations detected in our observational dataset, and 
thus place them below our strict 0.85 threshold. Also no- 
ticeable in the middle panels of Figure |4] is the reduced 
spread of phase angles for the simulated occurrence plots. 
To investigate whether this is a direct consequence of 
the noise-free time series, noise consistent with the ob- 
servational dataset was added, and the occurrence plots 
recreated. The acquired dataset has a combined photon 
and read noise of ~1% with respect to the observational 
count rates. The resulting image signal-to- noise ratio, 
SNR, can be calculated as SNR = /i/cr, where \i is the 
signal mean and a is the standard deviation of the noise 
(|Schroederll2000l ) . A signal mean of ^5000 counts, com- 
bined with a dark frame standard deviation count of ^19, 
provides an observational SNR w 260. The addition of 
1% noise to the simulated dataset did not modify the 
occurrence plots by a significant margin. However, when 
the noise level was increased to 5%, the peak occurrence 
phase angle and oscillation period remained consistent 
with the noise- free time series, but with a significantly 



larger spread of phase angles (lower panels of Fig. |4]). A 
noise level of 5% provides a simulated SNR w 150, still 
well above the Rose threshold for feature disambigua- 
tion ( Rosel [T948h . Thus, simulated features remain well 
defined following the addition of 5% noise, yet display a 
marked increase in the spread of recovered phase angles 
between the two filter bandpasses. Under these condi- 
tions, 80% of MBP structures display a negative phase 
shift, compared with 79% for non-MBP regions, indicat- 
ing a closer resemblance to the observational time series. 
Again, these results imply that, even in the simulated 
domain, there is a preference for 4170 A continuum os- 
cillations to lead their G-band counterparts. 

An important point to raise is how these oscillations 
are generated. Often, oscillations are produced in sim- 
ulated tim e series through the addition of an input 
driver (e.g. iBogdan et aTll2003l : ICally Goossensl 120081 
iFedun et al.ll2011bl ). However, in this instance, no spe- 
cific driver has been added to our time series. This 
implies that longitudinal oscillations readily arise as a 
natural consequence of convective motions, and their 
interaction with t he magnetic field in t he solar pho- 
tosphe re (see e.g. Kitias hvili et al.l 120111 : iShelvag et al.l 
l2011af ). Such oscillati o ns ha ve been studied observation- 
ally by iGoode et al.l (| 19981 ). who suggested these sig- 
natures are the result of events occurring just beneath 
the photosphere, which subsequently convert into atmo- 
sphe ric longitudinal oscill ations. Theoretical approaches 
(e.g. Ska rtlien et al.l 12000) have indicated that transient 
acoustic waves are excited at the top of the convective 
zone, and immediately above convective overshoot re- 
gions. Thus, the detection of an abundant array of fre- 
quencies over a range of photospheric structures, includ- 
ing MBPs and intergranular lanes, coupled with the rapid 
development of these oscillatory phenomena, further in- 
dicates how the ubiquity of oscillations may be traced 
back to magnetoconvective processes occurring through- 
out the solar photosphere. 

3.2.1. Determination of formation heights 

The contribution functions for our G-band and 4170 A 
continuum filters are displayed in Figure [5j for both MBP 
and granular structures, in addition to an average of the 
entire fleld-of-view. Our numerical domain is positioned 
in such a way that the visible solar surface, roughly de- 
fined as the horizontal geometrical layer which is physi- 
cally close to the optical layer of radiation formation, is 
located approximately 600 km belo w the upper boundar y 
(i.e. at a domain height of 800 km; She lvag et al.ll2011bf ). 
Initially, a peak subtraction method was utilised to esti- 
mate the height separation between the two filter band- 
passes. This process involved establishing the height at 
which the 4170 A continuum exhibited its largest contri- 
bution, and subtracting this from the simultaneous (and 
co-spatial) peak contribution height for the G-band filter. 
This always resulted in positive distances, ranging from 
approximately km at the center of granules, to over 
200 km in some intergranular lanes and MBPs (upper- 
right panel of Fig. [6]). However, from Figure [5j it is 
clear that the contribution functions for both G-band 
and 4170 A continuum filters are considerably asymmet- 
ric. A prolonged tail, present in the G-band contribu- 
tion functions, suggests that more emphasis should be 
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Fig. 6. — The height separation between G-band and 4170 A continuum bandpasses using the center-of-gravity (upper left) and peak 
subtraction (upper right) methods. In both upper panels, the color scale represents the separation distance (in km) using each of the 
methods detailed in §[3] Simultaneous G-band and 4170 A continuum intensity images, normalised to their mean, are shown for reference 
in the lower left and right panels, respectively. All axes are distances along the solar surface, where 1 Mm = 1000 km ~ 1.38" . 



placed on these extended heights, in order to accurately 
represent potential height separations. Thus, a center- 
of-gravity (COG) method was chosen to provide a more 
robust calculation of the distance between the two filter 
bandpasses. For each bandpass, the COG was calculated 
from the corresponding contribution function (CF) using, 



COG 



Ei = o(CF,xi) 



(3) 



where i represents individual measurements along the 
simulation's vertical domain, up to a maximum height, 
z, corresponding to 1.4 Mm. The COG measurement 
for the 4170 A continuum bandpass was then subtracted 
from the simultaneous (and co-spatial) COG value de- 
rived for the G-band filter. Under this regime, all es- 
timated distances were again positive. However, now 
a minimum height separation of approximately 60 km 
exists, even for contribution functions acquired at the 
center of granular cells, where the previous method esti- 
mated a km separation. This is a direct result of the 
inclusion of the extended G-band contribution tail seen 
in Figure[5] The upper-left panel of Figure [6] displays a 2- 
dimensional representation of the estimated atmospheric 



height separation between G-band and 4170 A continuum 
filter bandpasses using the COG method. Interestingly, 
this approach does not provide the same maximum sepa- 
ration distance as found in the peak subtraction method. 
Using the COG approach, a maximum separation of ap- 
proximately 170 km was found for some intergranular 
lanes and MBP structures, compared to a distance of 
nearly 220 km using the peak separation method. Nev- 
ertheless, the two methods provide similar conclusions, 
whereby the G-band is formed above the 4170 A contin- 
uum, particularly in regions where strong magnetic field 
concentrations are expected, such as intergranular lanes 
and MBPs. In these locations, an average height sep- 
aration of ~55 km and ~83 km is found for the peak 
subtraction and COG methods, respectively. Further- 
more, it is interesting to see how, in the locations of 
MBPs, the peaks of the G-band and 4170 A continuum 
contribution functions lie below the layer corresponding 
to the visible sola r surface. Th i s is consistent with th e 
deta iled models of ISpruitl (|1976f) , ISchiissler etHI (|2003f ) , 
and iShelvag et al.l (|2004 120101 ) , who describe how the 
creation of MBPs is a result of the tenuous interplay be- 
tween hot photospheric plasma and evacuated magnetic 
flux tubes. 
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Fig. 7. — Power spectra of the G-band (upper left) and 4170 A continuum (upper right) observations as a function of frequency and 
spatial wavenumber. The color scale represents the power (in orders of magnitude) above the quiescent background. The lower panel shows 
the phase delay between the G-band and 4170 A continuum oscillations as a function of frequency and spatial wavenumber. The color scale 
is in degrees, where positive and negative values indicate that the G-band leads, and trails, the 4170 A continuum, respectively. Limiting 
curves for acoustic- and gravity- wave regimes are displayed as dashed lines, with an evanescent region sandwiched between the acoustic 
and gravity regimes. A similar figure, only for the simulated time series, is shown in Figure 



The formation height for each filter bandpass can also 
be estimated by applying these techniques to the aver- 
aged field-of-view contribution functions. Using the peak 
subtraction method, a formation height of ~85 km and 
~40 km above the visible surface (located at ~800 km 
in Fig. [5]) is found for G-band and 4170 A continuum fil- 
ters, respectively. When the COG method is utilised, 
the formation heights of G-band and 4170 A continuum 
bandpasses becomes modified slightly to ~100 km and 
~25 km, respectively. It must be noted that these are 
lower limits of the height of formation. Molecular and 
atomic line absorption profiles, which exist in the filter's 
spectral bandpasses, will contribute to the overall radi- 
ation absorption, thus shifting the total (continuum + 
line) contribution functions. As a result, this will push 
radiation formation ranges upwards towards the cooler 
regions of the solar photosphere. Consequently, while 
the G-band will continue to be formed above the 4170 A 
continuum, the exact formation heights will be based on 
the lower limits established above. This has important 
consequences for the phase angles determined between si- 
multaneous G-band and 4170 A continuum oscillations. 



Since the height separation between the 4170 A contin- 
uum and the G-band is always positive, a negative phase 
shift can be attributed to an upwardly propagating; wave, 
such as tho se studied numerically by Wedeme ver et all 
(20031. 12004] ). Conversely, a positive phase shift will sug- 
gest the presence of a downwardly propagating wave, 
with a potential explanation incorporating aspects of 
wave reflection. However, the discrepancy between the 
amount of reflected waves in the observational (~28% 
with positive phase angles) and simulated (~5% with 
positive phase angles) time series may be related to the 
upper boundary conditions applied to the numerical sim- 
ulations, in addition to noise aspects discussed above. 
The simulated upper boundary is at least partially re- 
flective, but rigidly corresponds to a constant height in 
the photosphere. Contrarily, in the real Sun, this re- 
flective layer can vary as a function of time and spatial 
position, potentially leading to a broader scatter of phase 
differences for the reflected wave. 

3.2.2. Oscillation frequency as a function of spatial size 

While it is clear from Figure |4] that a spread of oscil- 
lation frequencies is present in the data, it is important 
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Fig. 8. — Power spectra of the G-band (upper left) and 4170 A continuum (upper right) simulations as a function of frequency and spatial 
wavenumber. The color scale represents the power (in orders of magnitude) above the quiescent background. The lower panel shows the 
phase delay between the G-band and 4170 A continuum oscillations as a function of frequency and spatial wavenumber. The color scale 
is in degrees, where positive and negative values indicate that the G-band leads, and trails, the 4170 A continuum, respectively. Limiting 
curves for acoustic- and gravity- wave regimes are displayed as dashed lines, with an evanescent region sandwiched between the acoustic 
and gravity regimes. A similar figure, only for the observational time series, is shown in Figure 



to investigate how individual photospheric structures af- 
fect the dete cted frequency. F ollowing the procedures 
detailed in iKrijger et al.l ([20011 ) , power spectra for the 
G-band and 4170 A continuum observations and simula- 
tions are shown as a function of oscillation frequency and 
spatial wavenumber in the upper panels of Figures [7] & 
El respectively. In each spectra, the majority of oscilla- 
tory power is associated with low frequencies and large 
spatial scales, con sistent with solar p-mode oscillations 
([Duvall et al.lll988h . These lower- left portions of the ob- 
servational and simulated spectra reveal power approach- 
ing six and ten orders of magnitude above the quiescent 
background, respectively, highlighting the considerable 
power associated with these waves. Moving to smaller 
spatial scales (i.e. higher wavenumbers), it is clear that 
higher oscillation frequencies tend to be associated with 
these structures. This effect is particularly pronounced 
in the upper panels of Figure [7J where a diagonal band of 
high oscillatory power is seen to stretch diagonally away 
from the lower-left corner of the power spectra. This 
can be verified in Figure HJ where oscillations down to 
~110 s can be detected in small-scale MBP structures, 



compared to a ^160 s lower limit on typically larger, 
non-MBP regions. 

It is interesting to note that almost no oscillations with 
periodicities below 100 s are detected simultaneously in 
the G-band and 4170 A continuum datasets. This ap- 
plies to both the observational and simulated data, but 
may be related to the broad filt er widths used. Recently, 
Kneer & Bello Gonzalez (2011) detected upwardly prop- 
agating waves in the lower solar atmosphere with peri- 
odicities as short as 70 s. The authors utilised 2-D spec- 
tral imaging techniques to obtain a time series with an 
extremely narrow equivalent bandpass (~18 mA). Con- 
trarily, the contribution functions extracted from our nu- 
merical simulations extend over a range of atmospheric 
heights due to the broad nature of our G-band and 
4170 A continuum filters, potentially decreasing and/or 
smoothing wave signals with periods shorter than 100 s. 
This does not necessarily mean that waves with period- 
icities below 100 s do not occur individually in either 
the G-band or 4170 A continuum datasets. However, it 
does indicate that if these waves exist, we do not detect 
their propagation between atmospheric layers since no 
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Fig. 9. A 4170 A continuum intensity image (left), normalised to its mean, is displayed alongside the simultaneous pixelised sound 
speeds (right), where the color scale represents the derived sound speeds in kms -1 . The average sound speed for the entire simulated 
field-of-view is 10.1 kms" 1 . All axes are distances along the solar surface, where 1 Mm = 1000 km ~ 1.38" . 



coherent signals can be determined. The phase distribu- 
tion displayed in the lower panel of Figure [71 reveals a 
shift of ^—10° between G-band and 4170 A continuum 
oscillations, apparent down to a spatial scale of ~0.4" 
(wavenumber 15 arcsec -1 ). This spatial size is consis- 
tent with typical (roundish) MBP diameters, and forms 
an upper wavenumber limit beyond which the phase dis- 
tribution becomes dominated by noise. These measure- 
ments are in agreement with the occurrence plots pre- 
sented in Figure HJ and reiterate the presence of upward 
propagating wave phenomena. 

Examination of the lower panel of Figure [71 reveals a 
frequency "fork" for wavenumbers below 4 arcsec -1 (spa- 
tial size of approximately 1.6"). This "fork" constitutes 
a continual linear downward trend of upwardly propagat- 
ing waves (negative phases) for oscillation frequencies as 
a function of spatial wavenumber, in addition to a visible 
plateau at frequencies ranging from 4-7 mHz, which 
includes the 3 min («5.5 mHz) p-mode oscillation. A 
plateau suggests that the frequency of detected oscilla- 
tions is no longer dependent on the increasing size of the 
spatial dimension, at least up to the ^70" limit imposed 
by our field-of-view. As a result, upwardly propagating 
3 minute oscillations appear to be ubiquitous over gran- 
ular (and larger) spatial scales . This agrees with the the- 
oretical work of lErdelvi et al.l (j2007l ). who demonstrate 
how 3 minute oscillations are free to propagate upwards 
from the solar surface as a consequence of their frequency 
residing above the acoustic cut-off value. A 3 minute pe- 
riodicity overlaps wit h the solar p-mode sp ectrum, which 
has been shown by iMoretti et al.l (j2QQlh to exist as a 
result of coherent pressure drivers up to spatial scales 
consistent with supergranular cells (~50 ff ). Thus, our 
results indicate an abundance of upwardly propagating 
longitudinal waves which cover the majority of our field- 
of-view, and suggests these wave forms are a direct re- 
sult of a coherent (supergranular) p-mode driver. Inter- 
estingly, a similar phase plateau is not apparent in the 
simulated frequency-wavenumber diagram (lower panel 



of Fig. [8j). This may be a direct result of the reduced 
field-of-view of the simulated time series. With an over- 
all size of 16.5 x 16.5 arcsec 2 , compared to the obser- 
vational field-of-view of 69.3 x 69.1 arcsec 2 , information 
occurring on large spatial scales (i.e. low wavenumbers) 
will be lost. 

Evanescent waves are a phenomenon known to exist 
in the frequency-wavenumber diagram in locations of 
very low frequency and small wavenumbe r (Frazier 1968; 
Deubnerl Il975at iLites fc Chipman! 119791 : iDeubner et al I 
1990). These waves cover the global 5 minute oscil- 
lation (~3 mHz), and are not expected to show sig- 
nals of upward propagation. Positioned on either side 
of the evanescent region are locations in the frequency- 
waven umber diagram where acoustic and gravi ty waves 
thrive (jSouffrinl 119661 : IStixl [l97Ql : iDeubnerl fl975l ) . To in- 
vestigate which regime best matches our detected wave 
phenomena, theoretical curves incorporating acoustic, 
evanescent, and gravity wave regimes, were over-plotted 
in the phase diagrams of Figures & [8j In order to do 
this, the pressure, P, and density, p, were first extracted 
from our numerical simulations. Then, the sound speed, 
c s , can be determined for each grid cell in our numerical 
domain using, 



(4) 



where 7 is the ratio of specific heats, here taken as 5/3. 
The pixelised sound speeds were then convolved with our 
bandpass-specific contribution functions to determine an 
accurate representation of the sound speeds typically 
contained within our filtergrams. A two-dimensional 
representation of the derived sound speeds is shown in 
the right panel of Figure [9j with an spatial average of 
10.1 kms -1 . Using this value for the sound speed, and 
follow ing the steps outlined in iKneer fc Bello Gonzalezl 
(|2011f h we derived the limiting curves for acoustic- and 
gravity- wave regimes. These frequency-wavenumber con- 
tours are over-plotted in the phase diagrams of Figures [7] 
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&[8j with the locations of acoustic, evanescent, and grav- 
ity wave regimes annotated. Close inspection of the lower 
panel in Figure [71 reveals how there are no upwardly prop- 
agating waves (i.e. negative phase angles) in areas of the 
diagram with frequencies less than ~3 mHz (5 min) and 
wavenumbers less than ~4 arcsec -1 (1.6"), consistent 
with previous studies of evanescent waves. In this portion 
of the frequency- wavenumber diagram, only downwardly 
propagating waves appear to prevail. A similar trend 
is detected in the simulated time series (lower panel of 
Fig. [8]), where downwardly propagating waves dominate 
for wavenumbers less than ~4 arcsec -1 and frequencies 
less than ~2.5 mHz. Slight differences between the obser- 
vational and simulated frequency- wavenumber diagrams 
at these extreme limits of frequency and wavenumber 
may again be a result of the reduced field-of-view size 
and time series length for the simulated data. Never- 
theless, the observational and simulated phase distribu- 
tions (lower panels of Figs.[3&[8j) reveal identical trends, 
whereby propagating waves, best described as acoustic 
phenomena, with a similar phase angle increase in fre- 
quency with decreasing spatial size. 

3.2.3. Wave propagation velocity 

A phase difference of —8° (as found in our observa- 
tional time series), in conjunction with the most com- 
monly occurring periodicity of 390 s, results in a time 
lag of 8.7 s between 4170 A continuum and G-band at- 
mospheric heights. Utilising a traversed distance of 55 - 
83 km, estimated through a comparison of G-band and 
4170 A continuum MBP contribution functions in our ra- 
diative MHD simulations, an average propagation speed 
of 7.9 ± 1.6 kms -1 is found. This derived wave speed is 
close to typical photospheric sound speeds, and reveals 
how, at least in the photospheric regime, the detected 
oscillations are linear in nature. 

4. CONCLUDING REMARKS 

Here we present high-cadence observations and simu- 
lations of the solar photosphere. Each dataset demon- 
strates a wealth of oscillatory behaviour, with high con- 
centrations found in highly magnetic regions, such as 
magnetic bright points (MBPs), visible as periodic inten- 
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sity fluctuations with periods in the range 110 - 600 s. 
Using a high-resolution radiative magneto-hydrodynamic 
(MHD) code to simulate photospheric magnetoconvec- 
tion, we detect a wealth of upwardly propagating acous- 
tic waves over a range of frequencies and wavenumbers. 
Results of the MHD simulations are consistent with those 
found in our observational time series, confirming that 
current MHD simulations are able to accurately repli- 
cate propagating wave phenomena in the lower solar at- 
mosphere, without the need of an external driver. We 
suggest how the ubiquity of these oscillations may be 
traced back to magnetoconvective processes occurring 
in the upper layers of the Sun's convection zone. We 
are able to estimate the average height of formation of 
our G-band and 4170 A continuum filters by convolving 
the filter bandpasses with the radiative output of our 
simulation. This provides minimum G-band and 4170 A 
continuum formation heights of 100 km and 25 km, re- 
spectively. We find that longitudinal oscillations exhibit 
a dominant phase delay of —8° between G-band and 
4170 A continuum observations, suggesting the presence 
of upwardly propagating waves. Almost no propagating 
waves with periods less than 100 s are detected in either 
observational or simulated data. This does not neces- 
sarily mean that these high-frequency oscillations do not 
occur individually in either bandpass. However, it does 
indicate that if these waves exist, we do not detect their 
propagation between atmospheric layers since no coher- 
ent signals can be determined. More than 73% of MBPs 
(73% from observations, 96% from simulations) display 
upwardly propagating wave phenomena, helping to ex- 
plain the ubiquitous nature of MHD waves in the outer 
solar atmosphere. 
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